// Updated on June 24, 2022 to correct a typo, after the first online publication (Jan 11, 2022) 

use ../processed/Census_aggr.dta, clear

// QOB Reg (Age 40-49 only; exclude about 200 obs with missing qob)
gen age_yob = (year-1)-birthyr
keep if age_yob>=40 & birthqtr<=4

egen cs = group(birthyr bpl)
recode birthyr (1910/1919=1) (1920/1929=2) (1930/1939=3), gen(cohort) 

local iv
forvalues j=1/3{
	forvalues k=1/3{
		gen Q`k'C`j' = (birthqtr==`k' & cohort==`j')
		local iv `iv' Q`k'C`j'
	}
}

local xlvls
forvalues j=1/18{
	gen x`j' = (yschl>=`j')
	local xlvls `xlvls' x`j'
}

matrix A = J(19,8,.)
mat colnames A = w_pop w_pop_se w_ols w_ols_se w_iv w_iv_se b_ols b_ols_se
local ctrl i.birthyr i.bpl
local op "[pw=wgt], vce(cluster cs)"
gen d0 = yschl==0
reg d0 [pw=wgt]
matrix A[1,1] =  _b[_cons]
matrix A[1,2] = _se[_cons]	
qui forvalues j=1/18{
	gen d`j' = yschl==`j'
	reg d`j' [pw=wgt]
	matrix A[`j'+1,1] =  _b[_cons]
	matrix A[`j'+1,2] = _se[_cons]	
	reg x`j' yschl `ctrl' `op'
	matrix A[`j'+1,3] =  _b[yschl]
	matrix A[`j'+1,4] = _se[yschl]		
	ivregress 2sls x`j' (yschl=`iv') `ctrl' `op'
	matrix A[`j'+1,5] =  _b[yschl]
	matrix A[`j'+1,6] = _se[yschl]	
	reg logwage yschl `ctrl' if yschl==`j' | yschl==`j'-1 `op'
	matrix A[`j'+1,7] =  _b[yschl]
	matrix A[`j'+1,8] = _se[yschl]		
}
preserve
clear
svmat A, names(col)
gen x = _n-1
order x, first
export delimited using ../result/wgt_x_qob.csv, replace
restore

rename bpl statefip
merge m:1 statefip using ../data/region_division_state.dta
drop if _merge==2
assert _merge==3
rename statefip bpl

local wlist
forvalues j=1/4{
	gen region_`j' = region==`j'
	local wlist `wlist' region_`j'
}
gen y60 = year==1960
gen y70 = year==1970
gen y80 = year==1980
gen a40 = (age_yob>=40 & age_yob<=44)
gen a45 = (age_yob>=45 & age_yob<=49)
local wlist `wlist' y60 y70 y80 a40 a45

reg yschl `iv' `ctrl' [pw=wgt]
predict xhat, xb // This line had a typo in the first version (resid was used instead of xb)
reg yschl `ctrl' [pw=wgt]
predict yschl_res, resid
reg xhat `ctrl' [pw=wgt]
predict xhat_res, resid

matrix B = J(9,8,.)
mat colnames B = w_pop w_pop_se w_ols w_ols_se w_iv w_iv_se b_ols b_ols_se
local j=0
local op "[pw=wgt], vce(cluster cs)"
qui foreach v in `wlist'{
	local j = `j'+1
	reg `v' [pw=wgt]
	matrix B[`j',1] =  _b[_cons]
	matrix B[`j',2] = _se[_cons]	
	gen w_`v' = `v'*yschl_res
	reg w_`v' yschl_res `op'
	matrix B[`j',3] =  _b[yschl_res]
	matrix B[`j',4] = _se[yschl_res]		
	ivregress 2sls w_`v' (yschl_res=xhat_res) `op'
	matrix B[`j',5] =  _b[yschl_res]
	matrix B[`j',6] = _se[yschl_res]		
	reg logwage yschl `ctrl' if `v'==1 `op'
	matrix B[`j',7] =  _b[yschl]
	matrix B[`j',8] = _se[yschl]		
}
preserve
clear
svmat B, names(col)
gen grp = "a"
order grp, first
local i=0
foreach v in `wlist'{
	local i=`i'+1
	replace grp = `"`v'"' in `i'
}
export delimited using ../result/wgt_w_qob_revised.csv, replace
restore
